#################### Replication - PSRM ####################

## This R Script replicates the analysis for Figures 1-2, Appendix Tables B1-B3,
# Appendix C,  Appendix D, and Appendix F (Except F2) in "Differential Efficacy of Survey Incentives 
# Across Contexts: Experimental Evidence from Australia, India, and the United States" 


rm(list=ls())


# Set working directory if needed 

# Create Figures and Tables folders
dir.create("Figures")
dir.create("Tables")

# Packages
library(readstata13)
library(magrittr)
library(dplyr)
library(ggplot2)
library(readxl)
library(tidyr)



# Load datasets
ind <- read.dta13("ind_incentives_psrm.dta")
us <- read.dta13("us_incentives_psrm.dta")
aus <- read.dta13("aus_incentives_psrm.dta")



####0. For text in section 3.1 ####

# 0a. India - create variable for percent complete 
# NOTE: complete (=1 if they finished the survey; =0 if they did not finish the survey BUT saw the incentives treatment)

options(digits=4)

#sawtreatment
sum(ind$sawtreatment, na.rm = TRUE)/length(ind$sawtreatment)*100 # (note: denominator=number of observations)

#completed survey
sum(ind$complete, na.rm=TRUE)/length(ind$sawtreatment)*100 #(note: denominator=number of observations)


# 0b. U.S. 2015 Study

#Response rate
sum(us$response, na.rm = TRUE)/length(us$response)*100

#completion rate
sum(us$complete, na.rm = TRUE)/length(us$complete)*100


# 0c. Australia Study

#response rate
sum(aus$response, na.rm=TRUE)/length(aus$response)*100

#completion rate
sum(aus$complete, na.rm = TRUE)/length(aus$complete)*100




####1. For Figure 1: Response Rates, All Countries ####

## 1a. India Study

# Response Rates by treatment - Approximate ITT


treat1 <- ind$answeredone[ind$treatmentnew==1]
treat2 <- ind$answeredone[ind$treatmentnew==2]
treat3 <- ind$answeredone[ind$treatmentnew==3]
treat4 <- ind$answeredone[ind$treatmentnew==4]
treat5 <- ind$answeredone[ind$treatmentnew==5]

mean_r <- rbind(sum(treat1, na.rm = T)/(14586/5), sum(treat2, na.rm = T)/(14586/5),
                       sum(treat3, na.rm = T)/(14586/5),sum(treat4, na.rm = T)/(14586/5),
                       sum(treat5, na.rm = T)/(14586/5))

error <- c(qt(0.95,sum(!is.na(treat1))-1), qt(0.95,sum(!is.na(treat2))-1),qt(0.95,sum(!is.na(treat3))-1),
           qt(0.95,sum(!is.na(treat4))-1),qt(0.95,sum(!is.na(treat5))-1))
  
se_r <- rbind((sd(treat1, na.rm = T)/sqrt(2917.2)), (sd(treat2, na.rm = T)/sqrt(2917.2)), (sd(treat3,na.rm = T)/sqrt(2917.2)),
                     (sd(treat4,na.rm = T)/sqrt(2917.2)),(sd(treat5,na.rm = T)/sqrt(2917.2)))

high_r <- mean_r + error*se_r
low_r <- mean_r - error*se_r
lgth_r <- error*se_r

india_approx_res <- cbind.data.frame(mean_r, high_r, low_r, lgth_r)
india_approx_res$treat <- c(1,2,3,4,5)



## 1b. U.S. 2015 Study


# Aggregate by treatment 

response <- us %>% 
  group_by(treatment) %>% 
  summarise(mean_r = mean(response), sd_r = sd(response), n_r=sum(!is.na(treatment)))

response <- subset(response[1:6,])
response$se <- (response$sd_r/sqrt(response$n_r))
response$error <- qt(0.95,df=response$n_r-1)

response$high <- response$mean_r + response$error*response$se
response$low <- response$mean_r - response$error*response$se
response$lgth <- response$error*response$se

us_res <- response

## 1c. Australia Study 


response <- aus %>% 
  group_by(treatment) %>% 
  summarise (sum_r = sum(response, na.rm=T), mean_r = mean(response, na.rm=T), sd_r = sd(response, na.rm = T), n_r=sum(!is.na(treatment)))

response$se <- (response$sd_r/sqrt(response$n_r))
response$error <- qt(0.95,df=response$n_r-1)

response$high <- response$mean_r + response$error*response$se
response$low <- response$mean_r - response$error*response$se
response$lgth <- response$error*response$se

aus_res <- response


#1d. Combine response rate data from three countries for figure


#Prepare to combine dataset 
#Relabel India's incentive numbers so they're consistent with other countries
ind_approx_res2 <- india_approx_res
ind_approx_res2$treatment[india_approx_res$treat==5] <- 6
ind_approx_res2$treatment[india_approx_res$treat==4] <- 5
ind_approx_res2$treatment[india_approx_res$treat==3] <- 4
ind_approx_res2$treatment[india_approx_res$treat==1] <- 1
ind_approx_res2$treatment[india_approx_res$treat==2] <- 2
colnames(ind_approx_res2)[2] <- "high"
colnames(ind_approx_res2)[3] <- "low"
colnames(ind_approx_res2)[4] <- "lgth"
ind_approx_res2$country <- rep("India",5)
ind_approx_res2 <- ind_approx_res2[,c(1:4,6,7)]


us_response2 <- us_res
us_response2 <- us_response2[,c(2,7,8,9,1)]
us_response2$country <- rep("US",6)


aus_response2 <- aus_res
aus_response2 <- aus_response2[,c(3,8,9,10,1)]
aus_response2$country <- rep("Australia",6)



#combine datasets from all three countries
response_all <- rbind(aus_response2, ind_approx_res2, us_response2)

labels <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery","Big Lottery","Big+Small Lottery")

# Plot Response Rate Figure for all countries 
ggplot(data=response_all, aes(x=treatment, y=mean_r, fill=country)) +
  geom_bar(stat="identity", position="dodge") + geom_errorbar(aes(ymin=low, ymax=high), width=.2,
                                                              position=position_dodge(.9)) + theme_light() +
  theme(plot.title = element_text(hjust = 0.5), axis.text.x  = element_text(angle=0, vjust=0.5, size=10)) + 
  labs(title="Response Rate - All Countries", 
       x="Treatment", y = "Response Rate") + scale_x_discrete(limits=labels) + scale_fill_grey() + 
  theme_minimal() + guides(fill=guide_legend(title="Country"))

ggsave("Figures/Response Rate_AllCountries_Main.jpeg", width = 8, height = 8)





#### 2. For Figure 2: Completion Rates, All Countries ####

## 2a. India Study

treat1 <- ind$complete[ind$treatmentnew==1]
treat2 <- ind$complete[ind$treatmentnew==2]
treat3 <- ind$complete[ind$treatmentnew==3]
treat4 <- ind$complete[ind$treatmentnew==4]
treat5 <- ind$complete[ind$treatmentnew==5]


mean_comp <- rbind(sum(treat1, na.rm = T)/(14586/5), sum(treat2, na.rm = T)/(14586/5),
                   sum(treat3, na.rm = T)/(14586/5),sum(treat4, na.rm = T)/(14586/5),
                   sum(treat5, na.rm = T)/(14586/5))

error <- c(qt(0.95,sum(!is.na(treat1))-1), qt(0.95,sum(!is.na(treat2))-1),qt(0.95,sum(!is.na(treat3))-1),
           qt(0.95,sum(!is.na(treat4))-1),qt(0.95,sum(!is.na(treat5))-1))

se_comp <- rbind((sd(treat1, na.rm = T)/sqrt(2917.2)), (sd(treat2, na.rm = T)/sqrt(2917.2)), (sd(treat3,na.rm = T)/sqrt(2917.2)),
                 (sd(treat4,na.rm = T)/sqrt(2917.2)),(sd(treat5,na.rm = T)/sqrt(2917.2)))

high_comp <- mean_comp + error*se_comp
low_comp <- mean_comp - error*se_comp
lgth_comp <- error*se_comp

india_approx_comp <- cbind.data.frame(mean_comp, high_comp, low_comp, lgth_comp)
india_approx_comp$treat <- c(1,2,3,4,5)

# 2b. U.S. 2015 study


complete <- us %>% 
  group_by(treatment) %>% 
  summarise (mean_comp = mean(complete), sd_comp = sd(complete), n_comp=sum(!is.na(complete)))

complete <- subset(complete[1:6,])
complete$se <- (complete$sd_comp/sqrt(complete$n_comp))
complete$error <- qt(0.95,df=complete$n_comp-1)

complete$high_c <- complete$mean_comp + complete$error*complete$se
complete$low_c <- complete$mean_comp - complete$error*complete$se
complete$lgth_c <- complete$error*complete$se
us_complete <- complete

# 2c. Australia Study

complete <- aus %>% 
  group_by(treatment) %>% 
  summarise (sum_comp = sum(complete, na.rm=T), mean_comp = mean(complete, na.rm = T), sd_comp = sd(complete, na.rm = T), n_comp=sum(!is.na(complete)))

complete$se <- (complete$sd_comp/sqrt(complete$n_comp))
complete$error <- qt(0.95,df=complete$n_comp-1)

complete$high_c <- complete$mean_comp + complete$error*complete$se
complete$low_c <- complete$mean_comp - complete$error*complete$se
complete$lgth_c <- complete$error*complete$se

aus_complete <- complete


#2d. Combine completion rate data from three countries for figure

#Prepare to combine dataset 
#Relabel India's incentive numbers so they're consistent with other countries
ind_approx_comp2 <- india_approx_comp
ind_approx_comp2$treatment[india_approx_comp$treat==5] <- 6
ind_approx_comp2$treatment[india_approx_comp$treat==4] <- 5
ind_approx_comp2$treatment[india_approx_comp$treat==3] <- 4
ind_approx_comp2$treatment[india_approx_comp$treat==1] <- 1
ind_approx_comp2$treatment[india_approx_comp$treat==2] <- 2
colnames(ind_approx_comp2)[2] <- "high_c"
colnames(ind_approx_comp2)[3] <- "low_c"
colnames(ind_approx_comp2)[4] <- "lgth_c"
ind_approx_comp2$country <- rep("India",5)
ind_approx_comp2 <- ind_approx_comp2[,c(1:4,6,7)]


us_complete2 <- us_complete
us_complete2 <- us_complete2[,c(2,7,8,9,1)]
us_complete2$country <- rep("US",6)


aus_complete2 <- aus_complete
aus_complete2 <- aus_complete2[,c(3,8,9,10,1)]
aus_complete2$country <- rep("Australia",6)

complete_all <- rbind(aus_complete2, ind_approx_comp2, us_complete2)

# Plot Completion Rate Figure for all countries
ggplot(data=complete_all, aes(x=treatment, y=mean_comp, fill=country)) +
  geom_bar(stat="identity", position="dodge") + geom_errorbar(aes(ymin=low_c, ymax=high_c), width=.2,
                                                              position=position_dodge(.9)) + theme_light() +
  theme(plot.title = element_text(hjust = 0.5), axis.text.x  = element_text(angle=0, vjust=0.5, size=10)) + 
  labs(title="Completion Rate - All Countries", 
       x="Treatment", y = "Completion Rate") + scale_x_discrete(limits=labels) + scale_fill_grey() + 
  theme_minimal() + guides(fill=guide_legend(title="Country"))

ggsave("Figures/Completion Rate_AllCountries_Main.jpeg", width = 8, height = 8)


#### APPENDIX B ####
#### Table B1 - Descriptive Statistics At Assignment: Teach For Australia ####


mean.pop <- aus %>%
  summarise(across(c(admit, selection_day, female, european, age_proxy), list(mean = ~ mean(.x, na.rm = TRUE))))

mean.samp <- aus %>%
  filter(response==1) %>%
  summarise(across(c(admit, selection_day, female, european, age_proxy), list(mean = ~ mean(.x, na.rm = TRUE))))

n.pop <- aus %>%
  summarise(across(c(admit, selection_day, female, european, age_proxy), list(~ sum(!is.na(.x)))))

n.samp <- aus %>%
  filter(response==1) %>%
  summarise(across(c(admit, selection_day, female, european, age_proxy), list(~ sum(!is.na(.x)))))

# Transpose the summaries
mean.pop <- mean.pop %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "Mean:Population") 

mean.samp <- mean.samp %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "Mean:Sample") 

n.pop <- n.pop %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "N:Population") 

n.samp <- n.samp %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "N:Sample") 

combined_aus<- bind_cols(mean.pop, mean.samp, n.pop, n.samp)

tableb1 <- xtable::xtable(combined_aus, digits=3)
print(tableb1, file="Tables/tableb1.txt")


#### Table B2 - Descriptive Statistics At Assignment: Teach For India ####

mean.pop <- ind %>%
  summarise(across(c(cohort_master, accepted_master, female_master), list(mean = ~ mean(.x, na.rm = TRUE))))

mean.samp <- ind %>%
  filter(sawtreatment==1) %>%
  summarise(across(c(cohort_master, accepted_master, female_master), list(mean = ~ mean(.x, na.rm = TRUE))))


n.pop <- ind %>%
  summarise(across(c(cohort_master, accepted_master, female_master), list(~ sum(!is.na(.x)))))

n.samp <- ind %>%
  filter(sawtreatment==1) %>%
  summarise(across(c(cohort_master, accepted_master, female_master), list(~ sum(!is.na(.x)))))

# Transpose the summaries
mean.pop <- mean.pop %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "Mean:Population") 

mean.samp <- mean.samp %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "Mean:Sample") 

n.pop <- n.pop %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "N:Population") 

n.samp <- n.samp %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "N:Sample") 

combined_ind<- bind_cols(mean.pop, mean.samp, n.pop, n.samp)

tableb2 <- xtable::xtable(combined_ind, digits=3)
print(tableb2, file="Tables/tableb2.txt")


#### Table B3 - Descriptive Statistics At Assignment: Teach For America ####

mean.pop <- us %>%
  summarise(across(c(admit, matriculated4, gpa2, schoolselective, female3, white2, pell, parentalcollege2), list(mean = ~ mean(.x, na.rm = TRUE))))

mean.samp <- us %>%
  filter(response==1) %>%
  summarise(across(c(admit, matriculated4, gpa2, schoolselective, female3, white2, pell, parentalcollege2), list(mean = ~ mean(.x, na.rm = TRUE))))

n.pop <- us %>%
  summarise(across(c(admit, matriculated4, gpa2, schoolselective, female3, white2, pell, parentalcollege2), list(~ sum(!is.na(.x)))))

n.samp <- us %>%
  filter(response==1) %>%
  summarise(across(c(admit, matriculated4, gpa2, schoolselective, female3, white2, pell, parentalcollege2), list(~ sum(!is.na(.x)))))

# Transpose the summaries
mean.pop <- mean.pop %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "Mean:Population") 

mean.samp <- mean.samp %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "Mean:Sample") 

n.pop <- n.pop %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "N:Population") 

n.samp <- n.samp %>%
  pivot_longer(everything(), names_to = "Statistic", values_to = "N:Sample") 

combined_us <- bind_cols(mean.pop, mean.samp, n.pop, n.samp)

tableb3 <- xtable::xtable(combined_us, digits=3)
print(tableb3, file="Tables/tableb3.txt")



#### C1.Table C1 - T-Test, Approximate Intent to Treat Response Rate, India ####


# Response Rates by treatment - Approximate ITT

treat1 <- ind$answeredone[ind$treatmentnew==1] # Control
treat2 <- ind$answeredone[ind$treatmentnew==2] # $5 Charity
treat3 <- ind$answeredone[ind$treatmentnew==3] # Small Lottery
treat4 <- ind$answeredone[ind$treatmentnew==4] # Big Lottery
treat5 <- ind$answeredone[ind$treatmentnew==5] # Mixed Lottery

mean_r <- rbind(sum(treat1, na.rm = T)/(14586/5), sum(treat2, na.rm = T)/(14586/5),
                sum(treat3, na.rm = T)/(14586/5),sum(treat4, na.rm = T)/(14586/5),
                sum(treat5, na.rm = T)/(14586/5))

se_r <- rbind((sd(treat1, na.rm = T)/sqrt(2917.2)), (sd(treat2, na.rm = T)/sqrt(2917.2)), (sd(treat3,na.rm = T)/sqrt(2917.2)),
              (sd(treat4,na.rm = T)/sqrt(2917.2)),(sd(treat5,na.rm = T)/sqrt(2917.2)))

approx_res <- cbind.data.frame(mean_r, se_r)
approx_res$treat <- c(1,2,3,4,5)

# T Test Function

t_test <- function(m1, m2, s1, s2, df)
{
  t <- (m1-m2)/(sqrt(s1^2 + s2^2)) 
  dat <- c(m1, m2, t, 2*pt(-abs(t),df))
  names(dat) <- c("Mean 1", "Mean 2", "t", "p-value")
  return(dat)
}

cont_5 <- t_test(approx_res$mean_r[1], approx_res$mean_r[2], approx_res$se_r[1], approx_res$se_r[2],
                 14586-1)

cont_small <- t_test(approx_res$mean_r[1], approx_res$mean_r[3], approx_res$se_r[1], approx_res$se_r[3],
                     14586-1)

cont_big <- t_test(approx_res$mean_r[1], approx_res$mean_r[4], approx_res$se_r[1], approx_res$se_r[4],
                   14586-1)

cont_mixed <- t_test(approx_res$mean_r[1], approx_res$mean_r[5], approx_res$se_r[1], approx_res$se_r[5],
                     14586-1)

five_small <- t_test(approx_res$mean_r[2], approx_res$mean_r[3], approx_res$se_r[2], approx_res$se_r[3],
                     14586-1)

five_big <- t_test(approx_res$mean_r[2], approx_res$mean_r[4], approx_res$se_r[2], approx_res$se_r[4],
                   14586-1)

five_mixed <-  t_test(approx_res$mean_r[2], approx_res$mean_r[5], approx_res$se_r[2], approx_res$se_r[5],
                      14586-1)

small_big <-  t_test(approx_res$mean_r[3], approx_res$mean_r[4], approx_res$se_r[3], approx_res$se_r[4],
                     14586-1)

small_mixed <- t_test(approx_res$mean_r[3], approx_res$mean_r[5], approx_res$se_r[3], approx_res$se_r[5],
                      14586-1)


big_mixed <- t_test(approx_res$mean_r[4], approx_res$mean_r[5], approx_res$se_r[4], approx_res$se_r[5],
                    14586-1)


ttest_itt_resp <- rbind(cont_5, cont_small, cont_big, cont_mixed,
                        five_small, five_big, five_mixed,
                        small_big, small_mixed,
                        big_mixed)

tablec1 <- xtable::xtable(ttest_itt_resp, digits=3)
print(tablec1, file="Tables/tablec1.txt")


#### C2.Table C2 - T-Test, Approximate Intent to Treat Completion Rate, India ####


treat1 <- ind$complete[ind$treatmentnew==1] # Control
treat2 <- ind$complete[ind$treatmentnew==2] # $5 Charity
treat3 <- ind$complete[ind$treatmentnew==3] # Small Lottery
treat4 <- ind$complete[ind$treatmentnew==4] # Big Lottery
treat5 <- ind$complete[ind$treatmentnew==5] # Mixed Lottery

mean_comp <- rbind(sum(treat1, na.rm = T)/(14586/5), sum(treat2, na.rm = T)/(14586/5),
                   sum(treat3, na.rm = T)/(14586/5),sum(treat4, na.rm = T)/(14586/5),
                   sum(treat5, na.rm = T)/(14586/5))

error <- c(qt(0.95,sum(!is.na(treat1))-1), qt(0.95,sum(!is.na(treat2))-1),qt(0.95,sum(!is.na(treat3))-1),
           qt(0.95,sum(!is.na(treat4))-1),qt(0.95,sum(!is.na(treat5))-1))

se_comp <- rbind((sd(treat1, na.rm = T)/sqrt(2917.2)), (sd(treat2, na.rm = T)/sqrt(2917.2)), (sd(treat3,na.rm = T)/sqrt(2917.2)),
                 (sd(treat4,na.rm = T)/sqrt(2917.2)),(sd(treat5,na.rm = T)/sqrt(2917.2)))

approx_comp <- cbind.data.frame(mean_comp, se_comp)
approx_comp$treat <- c(1,2,3,4,5)


# T Test Function
t_test <- function(m1, m2, s1, s2, df)
{
  t <- (m1-m2)/(sqrt(s1^2 + s2^2)) 
  dat <- c(m1, m2, t, 2*pt(-abs(t),df))
  names(dat) <- c("Mean 1", "Mean 2", "t", "p-value")
  return(dat)
}

cont_5 <- t_test(approx_comp$mean_c[1], approx_comp$mean_c[2], approx_comp$se_c[1], approx_comp$se_c[2],
                 14586-1)

cont_small <- t_test(approx_comp$mean_c[1], approx_comp$mean_c[3], approx_comp$se_c[1], approx_comp$se_c[3],
                     14586-1)

cont_big <- t_test(approx_comp$mean_c[1], approx_comp$mean_c[4], approx_comp$se_c[1], approx_comp$se_c[4],
                   14586-1)

cont_mixed <- t_test(approx_comp$mean_c[1], approx_comp$mean_c[5], approx_comp$se_c[1], approx_comp$se_c[5],
                     14586-1)

five_small <- t_test(approx_comp$mean_c[2], approx_comp$mean_c[3], approx_comp$se_c[2], approx_comp$se_c[3],
                     14586-1)

five_big <- t_test(approx_comp$mean_c[2], approx_comp$mean_c[4], approx_comp$se_c[2], approx_comp$se_c[4],
                   14586-1)

five_mixed <-  t_test(approx_comp$mean_c[2], approx_comp$mean_c[5], approx_comp$se_c[2], approx_comp$se_c[5],
                      14586-1)

small_big <-  t_test(approx_comp$mean_c[3], approx_comp$mean_c[4], approx_comp$se_c[3], approx_comp$se_c[4],
                     14586-1)

small_mixed <- t_test(approx_comp$mean_c[3], approx_comp$mean_c[5], approx_comp$se_c[3], approx_comp$se_c[5],
                      14586-1)


big_mixed <- t_test(approx_comp$mean_c[4], approx_comp$mean_c[5], approx_comp$se_c[4], approx_comp$se_c[5],
                    14586-1)


ttest_itt_comp <- rbind(cont_5, cont_small, cont_big, cont_mixed,
                        five_small, five_big, five_mixed,
                        small_big, small_mixed,
                        big_mixed)

tablec2 <- xtable::xtable(ttest_itt_comp, digits=3)
print(tablec2, file="Tables/tablec2.txt")



#### AUSTRALIA 
#### Aus: Tables D1-D8, Australia Regression Tables ####

# Name treatments 

aus$treatment <- factor(aus$treatment, 
                        levels = c(1,2,3,4,5,6),
                        labels = c("control", "fivecharity", "tencharity", "hundred", "thousand", "mixed"))


#### Table D1: Australia Response Rates, Linear Probability Model for ITT #### 


library(tidyverse)
library(texreg)

formulae <- list(
  aus$response ~ as.factor(relevel(aus$treatment, "fivecharity")),
  aus$response ~ as.factor(relevel(aus$treatment, "tencharity")),
  aus$response ~ as.factor(relevel(aus$treatment, "hundred")),
  aus$response ~ as.factor(relevel(aus$treatment, "thousand")),
  aus$response ~ as.factor(relevel(aus$treatment, "mixed"))
)

response_mods <- map(formulae, ~ lm(.x, data = aus, model = FALSE))

library(stringr)
modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled1 <- texreg(response_mods, stars = 0.001,
       custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
       omit.coef = "\\(Intercept\\)",
       custom.coef.names = coefnames,
       caption.above = TRUE,
       caption = "Australia: Response Rates, Linear Probability Model with Controls for ITT",
       
       include.adjr = TRUE,
       include.rsquared = TRUE,
       include.rmse = TRUE,
       include.nobs = TRUE,
       include.fstat=TRUE,
       digits=3)

print(tabled1, file="Tables/tabled1.txt")


# Only keep variables that don't repeat in the final table


#### Table D2: Australia Response Rates, Linear Probability Model with Controls for ITT #### 


library(tidyverse)
library(texreg)

formulae <- list(
  aus$response ~ as.factor(relevel(aus$treatment, "fivecharity"))+aus$age_proxy,
  aus$response ~ as.factor(relevel(aus$treatment, "tencharity"))+aus$age_proxy,
  aus$response ~ as.factor(relevel(aus$treatment, "hundred"))+aus$age_proxy,
  aus$response ~ as.factor(relevel(aus$treatment, "thousand"))+aus$age_proxy,
  aus$response ~ as.factor(relevel(aus$treatment, "mixed"))+aus$age_proxy
)

response_mods <- map(formulae, ~ lm(.x, data = aus, model = FALSE))

library(stringr)
modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Age Proxy",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled2 <- texreg(response_mods, stars = c(0.01, 0.05, 0.1),
       custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
       omit.coef = "\\(Intercept\\)",
       custom.coef.names = coefnames,
       caption.above = TRUE,
       caption = "Australia: Response Rates, Linear Probability Model with Controls for ITT",
       
       include.adjr = TRUE,
       include.rsquared = TRUE,
       include.rmse = TRUE,
       include.nobs = TRUE,
       include.fstat=TRUE,
       digits=3)

print(tabled2, file="Tables/tabled2.txt")


# Only keep variables that don't repeat in the final table



#### Table D3: Australia Completion Rates, Linear Probability Model #### 


library(tidyverse)
library(texreg)

formulae <- list(
  aus$complete ~ as.factor(relevel(aus$treatment, "fivecharity")),
  aus$complete ~ as.factor(relevel(aus$treatment, "tencharity")),
  aus$complete ~ as.factor(relevel(aus$treatment, "hundred")),
  aus$complete ~ as.factor(relevel(aus$treatment, "thousand")),
  aus$complete ~ as.factor(relevel(aus$treatment, "mixed"))
)

response_mods <- map(formulae, ~ lm(.x, data = aus, model = FALSE))


library(stringr)
modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled3 <- texreg(response_mods, stars = c(0.01, 0.05, 0.1),
       custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
       omit.coef = "\\(Intercept\\)",
       custom.coef.names = coefnames,
       caption.above = TRUE,
       caption = "Australia: Completion Rates, Linear Probability Model for ITT",
       
       include.adjr = TRUE,
       include.rsquared = TRUE,
       include.rmse = TRUE,
       include.nobs = TRUE,
       include.fstat=TRUE,
       digits=3)

print(tabled3, file="Tables/tabled3.txt")


# Only keep variables that don't repeat in the final table



#### Table D4: Australia Completion Rates, Linear Probability Model with Controls for ITT #### 


library(tidyverse)
library(texreg)

formulae <- list(
  aus$complete ~ as.factor(relevel(aus$treatment, "fivecharity"))+aus$age_proxy,
  aus$complete ~ as.factor(relevel(aus$treatment, "tencharity"))+aus$age_proxy,
  aus$complete ~ as.factor(relevel(aus$treatment, "hundred"))+aus$age_proxy,
  aus$complete ~ as.factor(relevel(aus$treatment, "thousand"))+aus$age_proxy,
  aus$complete ~ as.factor(relevel(aus$treatment, "mixed"))+aus$age_proxy
)

response_mods <- map(formulae, ~ lm(.x, data = aus, model = FALSE))

library(stringr)
modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Age Proxy",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled4 <- texreg(response_mods, stars = c(0.01, 0.05, 0.1),
       custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
       omit.coef = "\\(Intercept\\)",
       custom.coef.names = coefnames,
       caption.above = TRUE,
       caption = "Australia: Completion Rates, Linear Probability Model with Controls for ITT",
       
       include.adjr = TRUE,
       include.rsquared = TRUE,
       include.rmse = TRUE,
       include.nobs = TRUE,
       include.fstat=TRUE,
       digits=3)

print(tabled4, file="Tables/tabled4.txt")


# Only keep variables that don't repeat in the final table



#### Table D5: Australia Response Rates, Marginal Effects for ITT #### 


library(tidyverse)
library(texreg)

formulae <- list(
  aus$response ~ as.factor(relevel(aus$treatment, "fivecharity")),
  aus$response ~ as.factor(relevel(aus$treatment, "tencharity")),
  aus$response ~ as.factor(relevel(aus$treatment, "hundred")),
  aus$response ~ as.factor(relevel(aus$treatment, "thousand")),
  aus$response ~ as.factor(relevel(aus$treatment, "mixed"))
)

response_mods <- map(formulae, ~ glm(.x, data = aus, model = FALSE))


library(stringr)
modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled5 <- texreg(response_mods, stars = c(0.01, 0.05, 0.1),
       custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
       omit.coef = "\\(Intercept\\)",
       custom.coef.names = coefnames,
       caption.above = TRUE,
       caption = "Australia: Response Rates, Marginal Effects Model for ITT",
       
       include.adjr = TRUE,
       include.rsquared = TRUE,
       include.rmse = TRUE,
       include.nobs = TRUE,
       include.fstat=TRUE,
       digits=3)

print(tabled5, file="Tables/tabled5.txt")

# Only keep variables that don't repeat in the final table



#### Table D6: Australia Response Rates, Marginal Effects for ITT #### 


library(tidyverse)
library(texreg)

formulae <- list(
  aus$response ~ as.factor(relevel(aus$treatment, "fivecharity"))+aus$age_proxy,
  aus$response ~ as.factor(relevel(aus$treatment, "tencharity"))+aus$age_proxy,
  aus$response ~ as.factor(relevel(aus$treatment, "hundred"))+aus$age_proxy,
  aus$response ~ as.factor(relevel(aus$treatment, "thousand"))+aus$age_proxy,
  aus$response ~ as.factor(relevel(aus$treatment, "mixed"))+aus$age_proxy
)

response_mods <- map(formulae, ~ glm(.x, data = aus, model = FALSE))

library(stringr)
modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Age Proxy",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled6 <- texreg(response_mods, stars = c(0.01, 0.05, 0.1),
       custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
       omit.coef = "\\(Intercept\\)",
       custom.coef.names = coefnames,
       caption.above = TRUE,
       caption = "Australia: Response Rates, Marginal Effects with Controls for ITT",
       
       include.adjr = TRUE,
       include.rsquared = TRUE,
       include.rmse = TRUE,
       include.nobs = TRUE,
       include.fstat=TRUE,
       digits=3)

print(tabled6, file="Tables/tabled6.txt")


# Only keep variables that don't repeat in the final table




#### Table D7: Australia Completion Rates, Marginal Effects for ITT #### 


library(tidyverse)
library(texreg)

formulae <- list(
  aus$complete ~ as.factor(relevel(aus$treatment, "fivecharity")),
  aus$complete ~ as.factor(relevel(aus$treatment, "tencharity")),
  aus$complete ~ as.factor(relevel(aus$treatment, "hundred")),
  aus$complete ~ as.factor(relevel(aus$treatment, "thousand")),
  aus$complete ~ as.factor(relevel(aus$treatment, "mixed"))
)

response_mods <- map(formulae, ~ glm(.x, data = aus, model = FALSE))


library(stringr)
modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery",
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
               "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled7 <- texreg(response_mods, stars = c(0.01, 0.05, 0.1),
       custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
       omit.coef = "\\(Intercept\\)",
       custom.coef.names = coefnames,
       caption.above = TRUE,
       caption = "Australia: Completion Rates, Marginal Effects for ITT",
       
       include.adjr = TRUE,
       include.rsquared = TRUE,
       include.rmse = TRUE,
       include.nobs = TRUE,
       include.fstat=TRUE,
       digits=3)

# Only keep variables that don't repeat in the final table
print(tabled7, file="Tables/tabled7.txt")



#### Table D8: Australia Completion Rates, Marginal Effects for ITT #### 

  
  library(tidyverse)
  library(texreg)
  
  formulae <- list(
    aus$complete ~ as.factor(relevel(aus$treatment, "fivecharity"))+aus$age_proxy,
    aus$complete ~ as.factor(relevel(aus$treatment, "tencharity"))+aus$age_proxy,
    aus$complete ~ as.factor(relevel(aus$treatment, "hundred"))+aus$age_proxy,
    aus$complete ~ as.factor(relevel(aus$treatment, "thousand"))+aus$age_proxy,
    aus$complete ~ as.factor(relevel(aus$treatment, "mixed"))+aus$age_proxy
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = aus, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Age Proxy",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled8 <-  texreg(response_mods, stars = c(0.01, 0.05, 0.1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "Australia: Completion Rates, Marginal Effects with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)

# Only keep variables that don't repeat in the final table

print(tabled8, file="Tables/tabled8.txt")

  
  
###### INDIA 
#### India: Tables D9-D12, India Regression Tables ####
  
  # Name treatments 
  
  ind$treatmentnew <- factor(ind$treatmentnew, 
                             levels = c(1,2,3,4,5),
                             labels = c("control", "charity", "hundred", "thousand", "mixed"))
  
  
#### Table D9: India Response Rates, Linear Probability Model for ToT #### 

  
  
  library(tidyverse)
  library(texreg)
  
  formulae <- list(
    answeredone ~ as.factor(relevel(treatmentnew, "charity")),
    answeredone ~ as.factor(relevel(treatmentnew, "hundred")),
    answeredone ~ as.factor(relevel(treatmentnew, "thousand")),
    answeredone ~ as.factor(relevel(treatmentnew, "mixed"))
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = ind, na.rm=T, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "Small Lottery", "Large Lottery",
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery")
tabled9 <-  texreg(response_mods, stars = c(0.01, 0.05, 0.1),
         custom.model.names = c("$5 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "India: Response Rates, Linear Probability Model for ToT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)

  # Only keep variables that don't repeat in the final table
  
print(tabled9, file="Tables/tabled9.txt")

  
    
#### Table D10: India Completion Rates, Linear Probability Model for ToT #### 
  
  
  #library(tidyverse)
  #library(texreg)
  
  formulae <- list(
    complete ~ as.factor(relevel(treatmentnew, "charity")),
    complete ~ as.factor(relevel(treatmentnew, "hundred")),
    complete ~ as.factor(relevel(treatmentnew, "thousand")),
    complete ~ as.factor(relevel(treatmentnew, "mixed"))
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = ind, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "Small Lottery", "Large Lottery",
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery")
 tabled10 <- texreg(response_mods, stars = c(0.01, 0.05, 0.1),
         custom.model.names = c("$5 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "India: Completion Rates, Linear Probability Model for ToT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
 print(tabled10, file="Tables/tabled10.txt")
 
  
  
#### Table D11: India Response Rates, Marginal Effects Model for ToT #### 
  
  #library(tidyverse)
  #library(texreg)
  
  formulae <- list(
    answeredone ~ as.factor(relevel(treatmentnew, "charity")),
    answeredone ~ as.factor(relevel(treatmentnew, "hundred")),
    answeredone ~ as.factor(relevel(treatmentnew, "thousand")),
    answeredone ~ as.factor(relevel(treatmentnew, "mixed"))
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = ind,  model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "Small Lottery", "Large Lottery",
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery")
tabled11 <-  texreg(response_mods, stars = c(0.01, 0.05, 0.1),
         custom.model.names = c("$5 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "India: Response Rates, Marginal Effects Model for ToT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tabled11, file="Tables/tabled11.txt")

  
  
#### Table D12: India Completion Rates, Marginal Effects for ToT #### 
  
  
  #library(tidyverse)
  #library(texreg)
  
  formulae <- list(
    complete ~ as.factor(relevel(treatmentnew, "charity")),
    complete ~ as.factor(relevel(treatmentnew, "hundred")),
    complete ~ as.factor(relevel(treatmentnew, "thousand")),
    complete ~ as.factor(relevel(treatmentnew, "mixed"))
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = ind, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "Small Lottery", "Large Lottery",
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery", 
                 "Control", "$5 Charity", "Small Lottery" ,"Large Lottery")
tabled12 <-  texreg(response_mods, stars = c(0.01, 0.05, 0.1),
         custom.model.names = c("$5 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "India: Completion Rates, Marginal Effects Model for ToT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table

print(tabled12, file="Tables/tabled12.txt")

  
  ###### USA 
#### US: Tables D13-D20, U.S. Regression Tables ####
  
  
  us$treatment <- factor(us$treatment, 
                         levels = c(1,2,3,4,5,6),
                         labels = c("control", "fivecharity","tencharity", "hundred", "thousand", "mixed"))
  
  ##Table D13: U.S. Response Rates, Linear Probability Model for IIT #### 
  
  
  
  formulae <- list(
    response ~ as.factor(relevel(treatment, "fivecharity")),
    response ~ as.factor(relevel(treatment, "tencharity")),
    response ~ as.factor(relevel(treatment, "hundred")),
    response ~ as.factor(relevel(treatment, "thousand")),
    response ~ as.factor(relevel(treatment, "mixed"))
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled13 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S.: Response Rates, Linear Probability Model for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
 
print(tabled13, file="Tables/tabled13.txt")


#### Table D14: U.S. Response Rates, Linear Probability Model with Controls for IIT #### 
  
  
  
  formulae <- list(
    response ~ as.factor(relevel(treatment, "fivecharity")) + matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    response ~ as.factor(relevel(treatment, "tencharity"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    response ~ as.factor(relevel(treatment, "hundred"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    response ~ as.factor(relevel(treatment, "thousand"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    response ~ as.factor(relevel(treatment, "mixed"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", 
                 "Matriculated", "GPA", "School Select.", "Female", "White","Pell Grant Recipient","Parents Attended College",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled14 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S.: Response Rates, Linear Probability Model with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tabled14, file="Tables/tabled14.txt")

  
  
#### Table D15: U.S. Completion Rates, Linear Probability Model for IIT #### 
  
  formulae <- list(
    complete ~ as.factor(relevel(treatment, "fivecharity")),
    complete ~ as.factor(relevel(treatment, "tencharity")),
    complete ~ as.factor(relevel(treatment, "hundred")),
    complete ~ as.factor(relevel(treatment, "thousand")),
    complete ~ as.factor(relevel(treatment, "mixed"))
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled15 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S.: Completion Rates, Linear Probability Model for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tabled15, file="Tables/tabled15.txt")


  ##Table D16: U.S. Completion Rates, Linear Probability Model with Controls for IIT #### 
  
  
  formulae <- list(
    complete ~ as.factor(relevel(treatment, "fivecharity")) + matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    complete ~ as.factor(relevel(treatment, "tencharity"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    complete ~ as.factor(relevel(treatment, "hundred"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    complete ~ as.factor(relevel(treatment, "thousand"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    complete ~ as.factor(relevel(treatment, "mixed"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", 
                 "Matriculated", "GPA", "School Select.", "Female", "White","Pell Grant Recipient","Parents Attended College",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled16 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S.: Completion Rates, Linear Probability Model with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table

print(tabled16, file="Tables/tabled16.txt")

  
#### Table D17: U.S. Response Rates, Marginal Effects for IIT #### 
  
  formulae <- list(
    response ~ as.factor(relevel(treatment, "fivecharity")),
    response ~ as.factor(relevel(treatment, "tencharity")),
    response ~ as.factor(relevel(treatment, "hundred")),
    response ~ as.factor(relevel(treatment, "thousand")),
    response ~ as.factor(relevel(treatment, "mixed"))
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled17 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S.: Response Rates, Marginal Effects for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table

print(tabled17, file="Tables/tabled17.txt")


#### Table D18: U.S. Response Rates, Marginal Effects with Controls for IIT #### 
  
  
  
  formulae <- list(
    response ~ as.factor(relevel(treatment, "fivecharity")) + matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    response ~ as.factor(relevel(treatment, "tencharity"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    response ~ as.factor(relevel(treatment, "hundred"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    response ~ as.factor(relevel(treatment, "thousand"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    response ~ as.factor(relevel(treatment, "mixed"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", 
                 "Matriculated", "GPA", "School Select.", "Female", "White","Pell Grant Recipient","Parents Attended College",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled18 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S.: Response Rates, Marginal Effects with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tabled18, file="Tables/tabled18.txt")

  
  
#### Table D19: U.S. Completion Rates, Marginal Effects Model for IIT #### 
  
  formulae <- list(
    complete ~ as.factor(relevel(treatment, "fivecharity")),
    complete ~ as.factor(relevel(treatment, "tencharity")),
    complete ~ as.factor(relevel(treatment, "hundred")),
    complete ~ as.factor(relevel(treatment, "thousand")),
    complete ~ as.factor(relevel(treatment, "mixed"))
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled19 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S.: Completion Rates, Marginal Effects for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tabled19, file="Tables/tabled19.txt")


#### Table D20: U.S. Completion Rates, Marginal Effects with Controls for IIT #### 
  
  
  formulae <- list(
    complete ~ as.factor(relevel(treatment, "fivecharity")) + matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    complete ~ as.factor(relevel(treatment, "tencharity"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    complete ~ as.factor(relevel(treatment, "hundred"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    complete ~ as.factor(relevel(treatment, "thousand"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2,
    complete ~ as.factor(relevel(treatment, "mixed"))+matriculated4+gpa2+schoolselective+female3+white2+pell+parentalcollege2
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", 
                 "Matriculated", "GPA", "School Select.", "Female", "White","Pell Grant Recipient","Parents Attended College",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$10 Charity", "Small Lottery" ,"Large Lottery")
tabled20 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$10 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S.: Completion Rates, Marginal Effects with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
print(tabled20, file="Tables/tabled20.txt")

  
########## F. U.S. 2024 Replication
#### F. U.S: Tables F1,F3-F13, Figures F1-F2, U.S. 2024 Replication ####
  
  rm(list = ls())
  
  
  library(readstata13)
  us <- read.dta13("us2024_incentive_psrm.dta")  

  
  
  #### Figure F1: U.S. 2024 Study: Response Rates ####
  
  response <- us %>% 
    group_by(treatment) %>% 
    summarise (mean_r = mean(response), sd_r = sd(response), n_r=sum(!is.na(treatment)))
  
  response <- subset(response[1:6,])
  response$se <- (response$sd_r/sqrt(response$n_r))
  response$error <- qt(0.95,df=response$n_r-1)
  
  response$high <- response$mean_r + response$error*response$se
  response$low <- response$mean_r - response$error*response$se
  response$lgth <- response$error*response$se
  
  xtable::xtable(response, caption="US 2024 Replication: Response Rates", digits=3)
  
  library(ggplot2)
  ggplot(data=response, aes(x=treatment, y=mean_r)) +
    geom_bar(stat="identity") + geom_errorbar(aes(ymin=low, ymax=high), width=.2,
                                              position=position_dodge(.9)) + theme_light() +
    theme(plot.title = element_text(hjust = 0.5), axis.text.x  = element_text(angle=0, vjust=0.5, size=10)) + 
    labs(title="Response Rate, U.S.", 
         x="Treatment", y = "Response Rate") + scale_x_discrete(limits=labels)
  
  ggsave("Figures/US2024-Response Rate.jpeg", width = 8, height = 8)
  
  
  
  #### Figure F2: U.S. 2024 Study: Completion Rates ####
  
  complete <- us %>% 
    group_by(treatment) %>% 
    summarise (mean_comp = mean(complete), sd_comp = sd(complete), n_comp=sum(!is.na(complete)))
  
  complete <- subset(complete[1:6,])
  complete$se <- (complete$sd_comp/sqrt(complete$n_comp))
  complete$error <- qt(0.95,df=complete$n_comp-1)
  
  complete$high_c <- complete$mean_comp + complete$error*complete$se
  complete$low_c <- complete$mean_comp - complete$error*complete$se
  complete$lgth_c <- complete$error*complete$se
  
  xtable::xtable(complete, caption="US 2024 Replication: Completion Rates", digits=3)

  labels <- c("Control","$5 Charity", "$20 Charity", "Small Lottery", "Big Lottery", "Big+Small Lottery")
  
  ggplot(data=complete, aes(x=treatment, y=mean_comp)) +
    geom_bar(stat="identity") + geom_errorbar(aes(ymin=low_c, ymax=high_c), width=.2,
                                              position=position_dodge(.9)) + theme_light() +
    theme(plot.title = element_text(hjust = 0.5), axis.text.x  = element_text(angle=0, vjust=0.5, size=10)) + 
    labs(title="Completion Rate, U.S.", 
         x="Treatment", y = "Completion Rate") + scale_x_discrete(limits=labels)
  
  ggsave("Figures/US2024-Completion Rate.jpeg", width = 8, height = 8)
  
  
  #### Table F1: U.S. 2024 Study: Descriptive Statistics
  
  mean.pop <- us %>%
    summarise(across(c(admitted, matriculated,cumulativegpa,selective,female,white,pellgrant,age), list(mean = ~ mean(.x, na.rm = TRUE))))
  
  mean.samp <- us %>%
    filter(response==1) %>%
    summarise(across(c(admitted, matriculated,cumulativegpa,selective,female,white,pellgrant,age), list(mean = ~ mean(.x, na.rm = TRUE))))
  
  n.pop <- us %>%
    summarise(across(c(admitted, matriculated,cumulativegpa,selective,female,white,pellgrant,age), list(~ sum(!is.na(.x)))))
  
  n.samp <- us %>%
    filter(response==1) %>%
    summarise(across(c(admitted, matriculated,cumulativegpa,selective,female,white,pellgrant,age), list(~ sum(!is.na(.x)))))
  
  # Transpose the summaries
  mean.pop <- mean.pop %>%
    pivot_longer(everything(), names_to = "Statistic", values_to = "Mean:Population") 
  
  mean.samp <- mean.samp %>%
    pivot_longer(everything(), names_to = "Statistic", values_to = "Mean:Sample") 
  
  n.pop <- n.pop %>%
    pivot_longer(everything(), names_to = "Statistic", values_to = "N:Population") 
  
  n.samp <- n.samp %>%
    pivot_longer(everything(), names_to = "Statistic", values_to = "N:Sample") 
  
  combined_us<- bind_cols(mean.pop, mean.samp, n.pop, n.samp)
  
  tablef1 <- xtable::xtable(combined_us, digits=5)
  print(tablef1, file="Tables/tablef1.txt")
  
  
  #### Table F3: U.S. 2024 Study: Response and Completion Rates ####
  
  
  us$treatment <- factor(us$treatment, 
                         levels = c(1,2,3,4,5,6),
                         labels = c("control", "fivecharity", "twentycharity", "hundred", "thousand", "mixed"))
  
  table <- us %>% 
    group_by(treatment) %>% 
    summarise (mean_r = mean(response), sd_r = sd(response), 
               mean_c=mean(complete), sd_c=sd(complete),N=sum(!is.na(response)))
  
tablef3 <-  xtable::xtable(table, caption="US 2024 Replication: Response and Completion Rates", digits=3)

print(tablef3, file="Tables/tablef3.txt")


  
  
  #### Table F4: U.S. 2024 Study: Response Rates, Linear Probability Model for ITT ####
  
  formulae <- list(
    response ~ as.factor(relevel(treatment, "fivecharity")) ,
    response ~ as.factor(relevel(treatment, "twentycharity")),
    response ~ as.factor(relevel(treatment, "hundred")),
    response ~ as.factor(relevel(treatment, "thousand")),
    response ~ as.factor(relevel(treatment, "mixed"))
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", 
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery")
tablef4 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S. 2024 Study: Response Rates, Linear Probability Model for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tablef4, file="Tables/tablef4.txt")

  
  #### Table F5: U.S. 2024 Study: Response Rates, Linear Probability Model with Controls for ITT ####
  
  formulae <- list(
    response ~ as.factor(relevel(treatment, "fivecharity")) + age + cumulativegpa,
    response ~ as.factor(relevel(treatment, "twentycharity"))+age + cumulativegpa,
    response ~ as.factor(relevel(treatment, "hundred"))+age + cumulativegpa,
    response ~ as.factor(relevel(treatment, "thousand"))+age + cumulativegpa,
    response ~ as.factor(relevel(treatment, "mixed"))+age + cumulativegpa
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", 
                 "Age", "GPA",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery")
tablef5 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S. 2024 Study: Response Rates, Linear Probability Model with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tablef5, file="Tables/tablef5.txt")

  
  #### Table F6: U.S. 2024 Study: Completion Rates, Linear Probability Model for ITT ####
  
  formulae <- list(
    complete ~ as.factor(relevel(treatment, "fivecharity")) ,
    complete ~ as.factor(relevel(treatment, "twentycharity")),
    complete ~ as.factor(relevel(treatment, "hundred")),
    complete ~ as.factor(relevel(treatment, "thousand")),
    complete ~ as.factor(relevel(treatment, "mixed"))
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", 
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery")
tablef6 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S. 2024 Study: Completion Rates, Linear Probability Model for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tablef6, file="Tables/tablef6.txt")

  
  #### Table F7: U.S. 2024 Study: Completion Rates, Linear Probability Model with Controls for ITT ####
  
  formulae <- list(
    complete ~ as.factor(relevel(treatment, "fivecharity")) + age + cumulativegpa,
    complete ~ as.factor(relevel(treatment, "twentycharity"))+age + cumulativegpa,
    complete ~ as.factor(relevel(treatment, "hundred"))+age + cumulativegpa,
    complete ~ as.factor(relevel(treatment, "thousand"))+age + cumulativegpa,
    complete ~ as.factor(relevel(treatment, "mixed"))+age + cumulativegpa
  )
  
  response_mods <- map(formulae, ~ lm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", 
                 "Age", "GPA",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery")
tablef7 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S. 2024 Study: Completion Rates, Linear Probability Model with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
print(tablef7, file="Tables/tablef7.txt")

  
  
  #### Table F8: U.S. 2024 Study: Response Rates, Marginal Effects for ITT ####
  
  formulae <- list(
    response ~ as.factor(relevel(treatment, "fivecharity")) ,
    response ~ as.factor(relevel(treatment, "twentycharity")),
    response ~ as.factor(relevel(treatment, "hundred")),
    response ~ as.factor(relevel(treatment, "thousand")),
    response ~ as.factor(relevel(treatment, "mixed"))
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", 
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery")
tablef8 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S. 2024 Study: Response Rates, Marginal Effects for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
print(tablef8, file="Tables/tablef8.txt")

  
  
  #### Table F9: U.S. 2024 Study: Response Rates, Marginal Effects with Controls for ITT ####
  
  formulae <- list(
    response ~ as.factor(relevel(treatment, "fivecharity")) + age + cumulativegpa,
    response ~ as.factor(relevel(treatment, "twentycharity"))+age + cumulativegpa,
    response ~ as.factor(relevel(treatment, "hundred"))+age + cumulativegpa,
    response ~ as.factor(relevel(treatment, "thousand"))+age + cumulativegpa,
    response ~ as.factor(relevel(treatment, "mixed"))+age + cumulativegpa
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", 
                 "Age", "GPA",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery")
tablef9 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S. 2024 Study: Response Rates, Marginal Effects with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
  
print(tablef9, file="Tables/tablef9.txt")

  
  #### Table F10: U.S. 2024 Study: Completion Rates, Marginal Effects for ITT ####
  
  formulae <- list(
    complete ~ as.factor(relevel(treatment, "fivecharity")) ,
    complete ~ as.factor(relevel(treatment, "twentycharity")),
    complete ~ as.factor(relevel(treatment, "hundred")),
    complete ~ as.factor(relevel(treatment, "thousand")),
    complete ~ as.factor(relevel(treatment, "mixed"))
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", 
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery")
tablef10 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S. 2024 Study: Completion Rates, Marginal Effects for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table

print(tablef10, file="Tables/tablef10.txt")

  
  
  #### Table F11: U.S. 2024 Study: Completion Rates, Marginal Effects with Controls for ITT ####
  
  formulae <- list(
    complete ~ as.factor(relevel(treatment, "fivecharity")) + age + cumulativegpa,
    complete ~ as.factor(relevel(treatment, "twentycharity"))+age + cumulativegpa,
    complete ~ as.factor(relevel(treatment, "hundred"))+age + cumulativegpa,
    complete ~ as.factor(relevel(treatment, "thousand"))+age + cumulativegpa,
    complete ~ as.factor(relevel(treatment, "mixed"))+age + cumulativegpa
  )
  
  response_mods <- map(formulae, ~ glm(.x, data = us, model = FALSE))
  
  library(stringr)
  modelnames <- c("$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery", "Mixed Lottery")
  coefnames <- c("Control", "$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", 
                 "Age", "GPA",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery",
                 "Control", "$5 Charity", "$20 Charity", "Small Lottery" ,"Large Lottery")
tablef11 <-  texreg(response_mods, stars = c(0.01, .05, .1),
         custom.model.names = c("$5 Charity", "$20 Charity", "Small Lottery", "Large Lottery", "Mixed Lottery"),
         omit.coef = "\\(Intercept\\)",
         custom.coef.names = coefnames,
         caption.above = TRUE,
         caption = "U.S. 2024 Study: Completion Rates, Marginal Effects with Controls for ITT",
         
         include.adjr = TRUE,
         include.rsquared = TRUE,
         include.rmse = TRUE,
         include.nobs = TRUE,
         include.fstat=TRUE,
         digits=3)
  
  # Only keep variables that don't repeat in the final table
print(tablef11, file="Tables/tablef11.txt")

  
### Attention Check Questions, Tables F12-F14: Data Prep 
  
  library(dplyr)
  library(magrittr)
  library(tidyr)
  
  
  # Restrict to completed surveys
  
  
  # 10705 completed the survey
  # 10696 answered question
  
  completesonly <- subset(us, complete == 1)
  
#### Table F12: U.S. 2024 Study: Incentive Check Amongst Survey Completers ####
  
  check <- prop.table(table(completesonly$incentive_check, completesonly$treatment), margin=2)
tablef12 <-  xtable::xtable(check, caption = "U.S. 2024 Study: Incentive Check Amongst Survey Completers")
print(tablef12, file="Tables/tablef12.txt")

#### Table F13:  U.S. 2024 Study: Incentive Check Amongst Survey Completers, Charity Condition ####  
  check2 <- prop.table(table(completesonly$incentive_donation, completesonly$treatment), margin=1)
tablef13 <-  xtable::xtable(check2, caption="U.S. 2024 Study: Incentive Check Amongst Survey Completers, Charity Condition")
print(tablef13, file="Tables/tablef13.txt")

#### Table F14: U.S. 2024 Study: Incentive Check Amongst Survey Completers, Lottery Condition ####
  check3 <- prop.table(table(completesonly$incentive_lottery, completesonly$treatment), margin=1)
tablef14 <-  xtable::xtable(check3, caption = "Incentive Check Amongst Survey Completers, Lottery Condition")
print(tablef14, file="Tables/tablef14.txt")

  
  